%% Select local grid
idtemp=str2double(getenv('SGE_TASK_ID'));
if idtemp<workers
   u0grid_local=u0grid(jobs*(idtemp-1)+1: idtemp*jobs,:);
   u1grid_local=u1grid(jobs*(idtemp-1)+1: idtemp*jobs,:);
   u2grid_local=u2grid(jobs*(idtemp-1)+1: idtemp*jobs,:);
   u3grid_local=u3grid(jobs*(idtemp-1)+1: idtemp*jobs,:);
   u4grid_local=u4grid(jobs*(idtemp-1)+1: idtemp*jobs,:);
else
   u0grid_local=u0grid(jobs*(idtemp-1)+1: end,:);
   u1grid_local=u1grid(jobs*(idtemp-1)+1: end,:);
   u2grid_local=u2grid(jobs*(idtemp-1)+1: end,:);
   u3grid_local=u3grid(jobs*(idtemp-1)+1: end,:);
   u4grid_local=u4grid(jobs*(idtemp-1)+1: end,:);
end
sg_local=size(u1grid_local,1);


addpath /share/apps/mosek-9.2/9.2/toolbox/r2015aom

IdSetW=cell(n_years,1);

for j=1:n_years
    
load(num2str(all_years(j)));  
u1grid_local=V_Gumbel(1,x)*ones(sg_local,1);
PYX_cond_x=PXY_cond(x,:);


%% Components common across parameter values 
[sgreduced, ...
id_Ugridreduced, s_id_gridreduced, sU_id_gridreduced,...
P_col,...
Tcoord_zeros,...
Ugridreduced, n_U_sets, indices_pairs_extended_original, Tcoord, C, C_sign,...
cr_obsequiv, cr_one, cr_zeros,fill_obsequiv, fill_one, fill_zeros, numeqconstraints, RHS_one, RHS_zeros, RHS_obsequiv, ...
equalorder, n_vertices,infty_pos,...
cr_symmetry, fill_symmetry, RHS_symmetry, S1, S2, S3, ind,...
sU_id_gridreduced_2]=useful_anyparam(u0grid_local, u1grid_local, u2grid_local, u3grid_local, u4grid_local,sg_local,random_indices_pairs_all_R1{end}, random_indices_pairs_all_R2{end}, PYX_cond_x) ;  
       
old_length=7;     
indices_pairs_extended_restricted=cell(size(m_restricted_all_R2,2),1);
indices_pairs_extended_restricted{end}=indices_pairs_extended_original;
for t=1:size(m_restricted_all_R2,2)-1
    [C1, C2] = myind2ind (random_indices_pairs_all_R1{t}, old_length^n_U_sets);
    ind1=[C1 C2];
    [C1, C2] = myind2ind (random_indices_pairs_all_R2{t}, sU_id_gridreduced_2);
    ind2=[C1 C2];
    indices_pairs_extended_restricted{t}=[ind1; ind2];
end

clear indices_pairs_extended_original



%% Do the job
exists=zeros(sgreduced,1);
for g=1:sgreduced
    check=0;
    t=0;
    cont=1;
    
    while  cont==1 && check~=2 && t<=size(m_restricted_all_R2,2)
    t=t+1;
    [result, ...
     T, Tcoord,n_U_sets, n_vertices, C,idU1, idU2, idU3, idU4, s1,s2,s3,C_sign, sU,Aeq,beq]=...
         OneSide(param_MOSEK, g,id_Ugridreduced, s_id_gridreduced, sU_id_gridreduced,...
                 P_col,...
                 Tcoord_zeros,...
                 Ugridreduced, n_U_sets, indices_pairs_extended_restricted{t}, Tcoord, C, C_sign,...
                 cr_obsequiv, cr_one, cr_zeros,fill_obsequiv, fill_one, fill_zeros, numeqconstraints, RHS_one, RHS_zeros,RHS_obsequiv,n_vertices,infty_pos,...
                 cr_symmetry, fill_symmetry, RHS_symmetry, S1, S2, S3,...
                 ind);
            
    if strcmp(result.sol.itr.prosta, 'PRIMAL_AND_DUAL_FEASIBLE') && t<size(m_restricted_all_R2,2) 
       check=check_original_problem(tol,Ueq, Uineq,result.sol.itr.xx, T, Tcoord, indices_pairs_extended_restricted{end},n_U_sets, n_vertices, C, ...
                                    idU1, idU2, idU3, idU4, s1,s2,s3,C_sign, sU,...
                                    Aeq,beq);
       cont=1;                                          
       if check==2
          exists(g)=1;
       end
    elseif strcmp(result.sol.itr.prosta, 'PRIMAL_AND_DUAL_FEASIBLE') && t==size(m_restricted_all_R2,2)
        cont=0;
        exists(g)=1;
    else
        cont=0;
    end
    end
end
%% Enlarge results to the original grid of parameter vectors
exists=exists(equalorder);
IdSetW{j} = [u0grid_local(logical(exists),:)  u1grid_local(logical(exists),:)  u2grid_local(logical(exists),:) u3grid_local(logical(exists),:) u4grid_local(logical(exists),:)];


end



filename = sprintf('IdSetW.%d.mat', ID);
save(filename,'IdSetW')  
